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Abstract 

We use classical molecular dynamics simulations to study the collapse of single (SWNT) and 
double- walled (DWNT) carbon nanotube bundles under hydrostatic pressure. The collapse pressure 
(p c ) varies as 1/-R 3 , where R is the SWNT radius or the DWNT effective radius. The bundles show 
~ 30 % hysteresis and the hexagonally close packed lattice is completely restored on decompression. 
The p c of DWNT is found to be close to the sum of its values for the inner and the outer tubes 
considered separately as SWNT, demonstrating that the inner tube supports the outer tube and 
that the effective bending stiffness of DWNT, Ddwnt~2Dswnt- We use an elastica formulation 
to derive the scaling and the collapse behavior of DWNT and multi- walled carbon nanotubes. 

PACS numbers: 81.07.De, 02.70.Ns, 62.20.Dc, 62.50.+p 
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I. INTRODUCTION 



Since their discovery, carbon nanotubes have been subject to intense theoretical and 
experimental investigations due to their fascinating structural, electronic, and mechanical 
properties 1 . Carbon nanotubes are showing great promise in such diverse fields as nano- 
electronics, actuators, sensors^, nanofluidics, hydrogen storage, and high-strength materials. 
The mechanical properties of carbon nanotubes depend on the number of coaxial graphitic 
rings that go into their making. Significant advances have been made in the understanding 
of single (SWNT) and multi-walled (MWNT) carbon nanotubes. Double-walled carbon nan- 
otubes (DWNT) have been observed and synthesized*^, more recently. Being the simplest 
of the MWNT, they are ideal systems to study the evolution of various properties from the 
single to the multi-walled regime. 

High-pressure Raman experiments on SWNT bundles^^ point to a structural phase 
transition at ~ 2 GPa. The current understanding is that the initially circular nanotube 
cross section is distorted to an oval shape under pressure. High pressure X-ray diffrac- 
tion studies also indicate a phase transition^ from the ambient triangular lattice symmetry, 
which reappears under decompression. Molecular dynamics simulations suggest that SWNT 
bundles*^ as well as isolated tube a 11 ^ 2 collapse under hydrostatic pressure and that the 
collapse pressure varies as an inverse power law of the tube radius. More recently, several 
worker a 13 i 14 i 15 i 16 i 17 have used Raman spectroscopy to study bundles of DWNT under hydro- 
static pressure. They conclude that the environment inside the outer tube is highly defect 
free and unperturbed, that the outer tube acts as a protective shield for the inner tube and 
that the inner tube provides structural support to the outer tube. 

In this paper, we describe a set of molecular dynamics simulations performed to investi- 
gate the behavior of DWNT under pressure, focusing on the response of the inner and the 
outer tubes. These results are contrasted with similar MD simulations on bundles of SWNT. 
Observed results are interpreted within the framework of the elastica theory. 

II. SIMULATION METHODOLOGY 

We have used DREIDING^, a standard generic macromolecular force field, in all our 
molecular dynamics (MD) simulations. Table [I] lists the force field parameters used to 
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calculate intra and inter molecular interactions. Elliott et al. 9 have successfully used this 
force field to study the collapse of SWNT bundles under hydrostatic pressure. Our simu- 
lations have been performed using ModulaSim 19 , a modular and general purpose molecular 
modeling package. The ensemble used was one of constant particle number, pressure, and 
temperature (NPT). The temperature (300 K) and the applied hydrostatic pressure were 
maintained using the Berendsen thermostat and barostat^. The simulation cell consisted 
of 16 independent SWNT or DWNT arranged in a hexagonally close packed 4x4 bundle, 
with periodic boundary conditions and pressure applied along all three mutually perpen- 
dicular directions. The tubes were ten unit cells long (2.3 nm). It has been found 9 that 
nine independent tubes, ten unit cells long, are sufficient to avoid finite size effects. The 
MD simulations were carried out on four SWNT, (5,5), (10,10), (15,15), and (20,20) and 
four DWNT, (5,5)@(10,10), (7,7)0(12,12), (10,10)@(15,15), and (15,15)@(20,20) bundles 
using the standard velocity Verlet algorithm to integrate the equations of motion. The gap 
between the inner and the outer tubes is ~ 3.4 A, close to the inter-layer gap in graphite. 
The bundles were initially equilibrated at atmospheric pressure and subsequently subjected 
to step-wise monotonically increasing hydrostatic pressure increments, allowing the unit cell 
volume to equilibrate for at least 10 ps at each step. The simulation time step was 1 fs. 
Information about the structural transition was obtained by measuring the unit cell volume 
after equilibration at each hydrostatic pressure step. 



III. RESULTS AND DISCUSSION 



All the SWNT and DWNT equilibrated at atmospheric pressure have nearly circular cross 
sections, as shown in Figs. 1(a) and (c) for a SWNT and a DWNT bundle, respectively. 
At atmospheric pressure, r min /r max > 0.93 where r min and r max are the smallest and the 
largest distances from the center to the circumference of the tube cross section. To study the 
structural transition, we plot the reduced volume (V/Vq), where Vq is the unit cell volume 
at atmospheric pressure, for the various SWNT bundles as a function of pressure as shown 
in Fig. 2 (a). It is clear that each of the four SWNT bundles undergoes a spontaneous 
structural transition at a critical pressure (p c ), which decreases as the radius of the tubes 
increases, in agreement with previously published results*^. Unless otherwise specified, 
p c refers to the structural change pressure on the loading curve. On plotting the reduced 
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volume versus pressure for the DWNT bundles, as shown in Fig. 12(b), we once again observe 
clear structural transitions at well-defined critical pressures. Up to the critical pressure, the 
tube cross sections remain nearly circular with slight deformations from the circular shape. 
When the applied hydrostatic pressure exceeds p c , the tube cross sections assume an elliptical 
shape. Further increase in pressure results in a dumbbell shape as shown in Figs. 1(b) and 
(d). The loading and unloading curves show a ~ 30% hysteresis in all the bundles studied. 
The hysteresis is calculated as 100% x [ p l ° adin 9 _ p unioadin g y p ioadin 9 _ The hexagonally close 

packed lattice is completely restored in all SWNT and DWNT bundles on decompression. 
A closer look at the critical pressures of the DWNT bundles in Fig. El (b) reveals several 
interesting features. First, we notice that the p c of a DWNT bundle is greater than the p c 
of an SWNT bundle of the outer tubes alone. For example, the p c of the (10,10)@(15,15) 
DWNT is 4.1 GPa, a value higher as compared to the p c of (15,15) SWNT (0.9 GPa). This 
shows that the inner tube supports the outer tube under hydrostatic pressure. Second, the 
p c of the DWNT bundle is even higher than the p c of an SWNT bundle of the inner tubes 
alone. For the (10,10)@(15,15) tubes, 4.1 GPa (p c of the DWNT bundle) is higher than 3.2 
GPa (p c of (10,10) SWNT bundle). Having demonstrated that the p c of a DWNT bundle 
is higher than the p c of both the inner and the outer tubes, we now ask whether one can 
predict the p c of a DWNT bundle with the knowledge of the p c of the inner and the outer 
tubes. From Fig. |21 we see that the p c of a DWNT bundle is close to the sum of the p c of 
the inner and the outer tubes. In our example, 4.1 GPa (p c of the (10,10)@(15,15) DWNT 
bundle) is equal to 3.2 GPa (p c of (10,10) SWNT bundle) plus 0.9 GPa (p c of (15,15) SWNT 
bundle). In section IV, we derive an analytical result that demonstrates this behavior. 

Fig. El (a) shows the collapse pressure of the SWNT bundles versus the tube radius 
along with a 1/R 3 fi t 11 ' 12 . If we now define an effective radius of a MWNT with n walls, as 
■j^— = ~ XT=i W> cona P se pressure of the DWNT bundles is found to follow a 1/ R? e j^ 
dependence as seen from Fig. 01 (b). Let D be the bending modulus of the graphene sheet 
so that the energy per unit surface area associated with curvature k is given b y 21 ' 22 ' 2 ? 124 
u e = ^k 2 . The value of Dswnt has been estimated from a plot of the single point energies 
per unit surface area of seven isolated SWNT as a function of 1/R 2 (Fig. 0] (a)), which 
gives the value of Dswnt to be 2.90 eV. The mean curvature for a bundle was calculated 
by averaging its local value at each atom of each tube (see Appendix). Using the values of 
curvature and D (from Fig. 0J), the values of the elastic energy per unit surface area for all 
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the bundles were calculated and are plotted as a function of pressure in Fig. 5. It is clear 
that the elastic energy, as expected, follows the structural transition shown in Fig. El The 
insets show that the relative increase in elastic energy during collapse increases linearly with 
radius for both SWNT and DWNT, which can be understood by the following argument. 
The elastic energy per unit length of a nanotube of radius R before collapse is {7& c = Dir/R. 
After collapse, the tube has flat regions, which have no elastic energy, and bent regions. 
The shape of the bent regions is invariant with respect to change in radius R. That is, any 
increase in R simply increases the flat regions. Let R' be an effective radius of the bent 
regions. Then, the elastic energy of the collapsed nanotube is U ac = Dtt/R'. Since R' is a 
fixed number that is independent of R (but less than R in magnitude), the relative increase 
in elastic energy on collapse is 

Uac - Ubc _ ( Dtt/R' - Dn/R \ _ R 

U bc V Dn/R J R> U 

This shows that the percent change in elastic energy on collapse depends linearly on the 
radius of the nanotubes. 



IV. SCALING OF THE RESPONSE WITH RADIUS, ADHESION ENERGY, 
AND BENDING ELASTICITY. 

The main results so far, 1/R? dependence of p c and that the p c for DWNT is a sum 
of the p c for packing of separate SWNT bundles, can be derived from a formulation based 
on elastica theor y 2 ^ 2 ^ . Assume that the response of nanotube bundles to external loading 
can be calculated by minimizing a properly defined potential energy functiona l 23124 . For a 
nanotube bundle of length L, assume that deformations are primarily two-dimensional in 
the cross section of the bundle. Let s denote length along a path that traverses the graphene 
sheets in this 2D cross section. Let S represent the complete path and Sp the part of this 
path on which we apply external pressure. Tubes interact with each other via an interatomic 
potential that has an attractive van der Waals component and a short-range repulsion. The 
latter component of the interaction effectively prohibits two surfaces from approaching each 
other too closely. Therefore, the surface S can be considered as consisting of two parts. 
Over the part of S where the short-range repulsion creates a flat interface, S a , we say they 
are in contact; over the remaining surface, S — S a , we say that they are not in contact. The 
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utility of this partitioning is that a material property, the work of adhesion of two graphene 
sheets, can be associated with S a . The potential energy per unit length of the bundle can 
be written a a 23 i 24 

Y = [ ^-k 2 ds - W a S a - [ u e v f w ds -I T-vds+ [ u g ds, (2) 

L JS 1 JS-Sa JS T JS 

where the first term represents energy due to bending of the graphene sheet to mean curva- 
ture k, the second and third terms capture adhesive van der Waals interactions, the fourth 
term represents the work of external forces with T, the external traction, and v, the displace- 
ment on the surface where tractions are applied, and the last term is the energy of formation 
of a flat sheet (u g is the energy of formation per unit surface area of a graphene sheet). The 
attractive van der Waals energy has been written in two parts. The first, W a S a , captures 
regions where graphene sheets are in contact; W a is the work of adhesion per unit area of 
bringing two flat nanotube walls from infinity to equilibrium separation. For DWNT and 
MWNT, there is a contribution to the work of adhesion due to interlayer contact. Because 
layers deform together, this contribution does not change with deformation. It will therefore 
vanish in a variation and for DWNT/MWNT S a can be identified as the area of contact 
between outermost layers. The term u^ w represents interactions outside the contacting 
regions. Once the tubes are in contact, the change in this term with further deformation 
can be neglecte d 23 ! 24 . The scaling of the solution can be extracted simply by a suitable 
normalization. Normalizing all length scales in Eq. |2] by the radius, and dropping terms 
that vanish in a variation, we obtain an expression for potential energy per unit length, v, 

Rv r 1 r2 = r r _ J _ w a R 2 r tr 3 a n pr? 

v = — = J -k ds - ab a + J a • vcis; a = — ; b = — — = fjn; (j = (3) 

where b is a dimensionless applied traction field, (3 is its (scalar) value for the case of 
fixed applied pressure, and n is the unit normal. Our dimensionless formulation implicitly 
assumes that no other length scale enters into the problem, for example, through boundary 
conditions. A possible exception is the interlayer spacing. Prior to collapse, it has been 
shown that the deformation of two nanotubes in contact is independent of this paramete r 23 ! 24 . 
In the collapsed state the interlayer spacing perturbs only slightly the solution obtained by 
neglecting i t 23 ! 24 . On this basis, we neglect its influence on our formulation; this assumption 
is justified by good agreement between the predicted scaling and simulation results. 

Therefore, within our assumptions, the deformation depends only on two dimensionless 
parameters, a and (3. If the dominant influence on deformation is external pressure, then 



events such as collapse or phase transitions will occur at critical values of j3, say at (3 C . This 
establishes the scaling of critical pressure to be P c = /3 C D/R 3 . In the absence of external 
pressure such events will occur at critical values of a, as already establishe d 23 ^ 4 . 

For commensurately packed MWNT, where difference in radius equals the equilibrium 
separation between graphene sheets, this argument can be extended to DWNT and MWNT. 
First, we recognize that the shape of any one shell in a deformed MWNT can be obtained 
from the shape of another shell simply by a change of scale. To build up a deformed 
MWNT packing, we therefore first start with a SWNT packing, say the shell with the largest 
diameter. Consider Fig. 1 (a), a SWNT bundle, in the case of deformations dominated by 
external pressure. At any stage of the deformation the solution A; is a function of a, (3. Denote 
by T the surface tractions needed to support the shape of this shell. Now make an identical 
copy of the deformed bundle and reduce its diameter by a change in scale to a value just 
small enough to fit inside the first shell. The new shape is also a solution if we scale all the 
tractions according to T = f3D/R 3 . Take the smaller bundle, separate the nanotubes, and 
insert each into its corresponding tube in the larger bundle. If we assume that the interface 
between the walls of a DWNT cannot carry any shear tractions, the resulting bundle is also 
a solution. We note that the tractions on the inner tube have to be provided by the outer 
tube. This leads to the conclusion that the net external tractions we need to apply to the 
outer shell is the sum of tractions needed to bring the two constituent SWNT's to similar 
shapes. This is easily generalized to a MWNT, establishing the fact that for deformation to 
the same normalized shape, the needed applied pressure is p = YH=iPi- AH the nanotubes 
will collapse simultaneously and so 

n 1 n 1 

Pc = Y.P™ = ^ nD ~ ^3 ( 4 ) 

i=l i % 

thus establishing the result that the collapse pressure for a MWNT packing is the sum of 
collapse pressures of the constituent SWNT packings, and providing the rationale for the 
effective radius defined earlier. 

As an independent test of this model, we have plotted the collapse pressure as a function 
of R e ff in Fig. 3 (b); it fits a l/-^e// relationship well. The fit yields a value for /3 c D e ff = 8.92 
eV. If we define the surface of DWNT to be a cylinder with radius R e ff, the energy per unit 
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surface area of DWNT is given by 

u e = _ 5 ( ^- + ^ - W«7r(i2i + R 2 ) + 2iru g {Ri + R 2 ] 

where 1/i? 2 = (l/2RiR e ff + l/2R 2 R e ff), Ri is the radius of the inner tube, and R 2 is the 
radius of the outer tube. In Eq. the first term corresponds to the elastic energies of the 
inner and outer tubes, the second term to the interaction energy of the two tubes and the 
third term to the energy of formation of the two tubes. Eq. El shows that the energy per 
unit surface area of DWNT scales with the inverse square of the length R. This quantity is 
readily computed for different DWNT and Fig. 4(b) plots the single point energy per unit 
surface area of five DWNT as a function of 1/R 2 . A fit using22*24 Ug = 0.765 and W a = 0.4 
yields a value of D e ff = 6.4 eV, which is close to twice Dswnt- Together with the fit to 
collapse pressure, the DWNT data yield a value of p^ WNT = 1.39, very close in value to 
the that obtained from SWNT simulations ((3^ WNT = 1.31). 



V. SUMMARY 



To summarize, we use classical MD simulations to show that DWNT bundles collapse 
at a critical pressure p c that, like in the case of SWNT, varies as l/R^f*, where R e ff is a 
suitably defined effective radius. We find that the SWNT and DWNT bundles show a ~ 
30% hysteresis and that the hexagonally close packed lattice is completely restored in all 
SWNT and DWNT bundles on decompression. Interestingly, we find that the p c of a DWNT 
bundle varies as the sum of the p c of the inner and the outer tubes considered separately as 
SWNT bundles (a result we derive analytically), demonstrating that the inner tube supports 
the outer tube and that D dwnt ~2D SW nt, where D is a bending stiffness. 
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Note. After the completion of our analysis, Ye et al.— published constant pressure MD 
simulations demonstrating a hydrostatic pressure-induced structural transition for isolated 
DWNT. The values of the critical pressures they obtain for isolated DWNT are 0.4 to 0.5 
times the values we find for the same diameter DWNT arranged in a bundle. 

APPENDIX: Calculation of mean curvatures for SWNT and DWNT bundles 

The following algorithm was used for the calculation of the mean curvature for a bundle 
of SWNT or DWNT. Each of the 16 tubes in the system is an armchair tube (n, n) with 
ten unit cells. It can be shown that the total number of atoms per tube is 40n. The local 
curvature is calculated at every atom that belongs to the middle eight unit cells (32n atoms). 
The atoms belonging to the unit cells at the ends of the tubes are not considered because 
these atoms do not have the sufficient number of neighbors required for our calculations (as 
will be clear later). 

Each of the 32n atoms is considered one at a time. For each atom, the coordinates of its 
three nearest neighbors and six next nearest neighbors are found using a search algorithm. 
The central atom's three nearest neighbors are used to define a plane passing through them 
and the normal to this plane is found. This is defined to be the new z-axis. The new x- and 
y-axes are suitably defined to be mutually perpendicular. 

A rotation matrix is now constructed using the components of the normal. The matrix 
is then used to transform the coordinates of the ten atoms (the central atom and its nine 
neighbors) to the new coordinate system. In the new coordinate system, a quadratic surface 
of the form z = g(x, y) is fit to the ten points as follows. The expanded form of the equation is 
given by z = ax 2 +by 2 +cxy+dx+ey+f. This equation can be written treating (a, b, c, d, e, /) 
as the unknowns and (x 2 , y 2 , xy, x, y, z) as the coefficients. The coordinates of the ten points 
give us ten equations in six unknowns. In matrix notation, we have [A]i 0X 6^4]6xi = [^]ioxi, 
where [A] is the matrix to be determined. The values of (a, b, c, d, e, /) are obtained by 
calculating [A] using the relation [A] = ( [N] T [N] ) ~ x [A] T [Z\ . 

The mean curvature of a surface, as defined above, is given by^ 




(6) 
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The value of H at the central atom is now calculated using the values of (a, b, c, d, e, /) and 
the coordinates of the central atom. 

This process is repeated for all the 32n atoms of the tube. The local curvature values 
at atoms of the other fifteen tubes in the bundle are similarly calculated to yield a total of 
512n values. The average curvature for the bundle is simply the mean of these 512n values. 
For bundles of DWNT, the same procedure is used, treating the inner and the outer tubes 
as separate SWNT and averaging over atoms in 32 tubes. 

This method gave good results for all tubes except the very small (5,5) tube, which cannot 
be well approximated by a smooth cylinder even at K. The calculated mean curvature, 
using the method described above, for an optimized (5,5) tube at K differs from the value 
of the curvature of a cylinder of the same radius (given by 1/ radius) by more than 5 %. 
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Ebond(R) = \K b {R - R ) 2 R 1.39 A K b 1050 (kcal/mol) /A 2 

E ang ie(0) = ±K (cos 6 - cos 6 ) 2 O 120° K e 100 (kcal/mol)/rad 2 

E t orsion{<t>) = \V{1 - cos[n(0 - O )]} <t>o 180° V 25.0 kcal/mol n 2 

#in«(*) = 3 (il^M™ 8 * " cos ^o) 2 #0 0° if* 40 (kcal/mol) /r ad 2 

E v dw(R) = Do {(^) 12 - 2 (I 1 ) 6 } go 3 - 8983 °j_ D o 0- 0951 kcal/mol 

TABLE I: Parameters for the C_R atom type (sp 2 hybridized carbon atom involved in resonance), 
in DREIDING 1 -, a standard generic macromolecular force field used in all our molecular dynamics 
simulations. 



12 




FIG. 1: The upper figures show a 4 x 4 bundle of (10,10) SWNT at (a) p = 1.0 Atm (before 
collapse) and (b) p = 6.0 GPa (after collapse). The bottom figures show a 4 x 4 bundle of 
(10,10)0(15,15) DWNT at (c) p = 1 Atm (before collapse) and (d) p = 10.2 Gpa (after collapse). 
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FIG. 2: Reduced volume (V/Vq) as a function of applied hydrostatic pressure for (a) SWNT and 
(b) DWNT bundles. The loading (solid symbols) and unloading (open symbols) curves clearly 
show hysteresis. 
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FIG. 3: Critical collapse pressure (p c ) as a function of (a) SWNT radius and (b) DWNT effective 
radius denned in the text. We estimate the values of f3 c for both SWNT and DWNT using Eq. |3J 
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FIG. 4: Single point energy per unit surface area at K as a function of 1/R 2 for seven isolated 
SWNT (a) and as a function of 1/R 2 for five isolated DWNT (b). The value of D, the bending 
stiffness, is obtained from the fit to Eq. [21 
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FIG. 5: Energy per unit surface function of applied hydrostatic pressure for (a) SWNT, 

and (b) DWNT bundles. Notice the correspondence with Fig. |2 The insets show the relative 
increase in elastic energy during collapse as a function of radius. 
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